Mon. Not. R. Astron. Soc. 000, 000-000 (0000) Printed 29 August 201 1 (MN KTeX style file v2.2) 



Portrait of Theobalda as a Young Asteroid Family 



Bojan Novakovic * 

Department of Astronomy, Faculty of Mathematics, University of Belgrade, Studentski trg 16, 11000 Belgrade, Serbia 



o 

(N 

o 



Oh; 

6 



(N 
> 

rn 

o 
o 



% 



29 August 2011 



ABSTRACT 

The (778) Theobalda asteroid family attracted little attention so far, but our study shows that 
it is important in several aspects. In this paper we investigate the origin and evolution of 
Theobalda family. Firstly, we identify the family as a statistically relevant group in the space 
of synthetic proper elements. Using the hierarchical clustering method and adopted cut-off 
velocity of dcutoff = 85 ms~^ we found that Theobalda family currently consists of 128 
members. This family is located in the outer belt, near proper semi-major axis Cp « 3. 175 au. 
This region is crossed by several three-body mean motion resonances which give rise to signif- 
icant chaotic zones. Consequently, the majority of family members reside on chaotic orbits. 
Using two independent methods, chaotic chronology and backward integration, we found 
Theobalda family to be only 6.9 ± 2.3 Myr old. We have also estimated, that the family was 
likely produced by the cratering impact on a parent body of diameter Dps ~ 78 ± 9 km. 

Key words: celestial mechanics, minor planets, asteroids, methods: numerical 



1 INTRODUCTION 



Asteroid families are believed to originate by a catastrophic disrup- 
tion of large asteroids. To identify an asteroid family, one looks 
. for clusters of asteroids i n the space of proper orbital elements 
' JMilani & Knezeviq|l990l . ll994}) : the proper semi-major axis (ap), 
, proper eccentricity (Cp) and proper inclination (Ip). The orbital el- 
' ements describe the size, shape and tilt of orbits. Proper orbital el- 
ements, being more constant over time than instantaneous ones, 
provide a dynamical criterion of whether or not a group of bodies 
' has a common ancestor. Up to now, ejecta from a few ten s of major 
I collisions have been discovered in the main belt (e.g. Zap pala et al.l 
! ll994lMothe-Diniz et al.ll2005l) . 

The size and velocity distributions of the family members 
provide important constraints for testing our understanding of the 
break-up process, but erosion and dynamical evolution of the or- 
bits over time can alter the original signature of the collision. It is 
nowadays well known that the kinematical structures of the asteroid 
families evolved over the time, with respect to the original post- 
impact situations, due to c haotic diffusion, gravitational and non- 
gravitational perturbations ( Milani & Farinella 19941: ISottke et al.l 
[2001; Carrubaetal. 2003; Dell'Oro et al. 2004). These mecha- 
nisms changed the original shapes of the families produced in col- 
lisions, and consequently complicated physical studies of high- 
velocity collisions. 

Unfortunately, most of the ob served astero id families are old 
enough (older than 100 Myr ( Nesvornv et al.l i2006 )) to be sub- 
stantially eroded and dispersed. On the other hand, young asteroid 
families (younger than 10 M yr) such as Karin, Veritas and lannini 
JNesvomv et alj2002ll2003h or even very young families (younger 
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than 1 Myr) such as Datura, Emi lkowalski, 1992YC2, and Lucas- 
cavin JNesvomv & Vokrouhlickv:.2006.) , suffer little erosion dur- 
ing the period of time after a breakup event. Thus, they provide 
a unique opportunity to study a collisional outcome almost unaf- 
fected by orbit evolution. 

In this paper we study Theobalda asteroid family. We present 
its basic properties including the identification of its membership, 
and the study of cumulative absolute magnitude distribution of the 
family members. Moreover, the diameter of the parent body has 
been estimated. We studied in detail the dynamical characteristics 
in the region occupied by the Theobalda asteroid family and ana- 
lyzed the role of the dynamics in shaping the family . 

As was noted by iNovakovic et al.l ( l2010ab this family is a 
very good candidate to estimate its age by the method of chaotic 
chronology (MCC). In order to apply MCC the family has to be 
located in the region of the main asteroid belt where diffusion 
takes place. Also, it is necessary that diffusion is fast enough 
to cause measurable effects, but slow enough so that most of 
the family members are still forming a robust family structure. 
As we show, it turns out that Theobalda family is an excellent 
case in this respect. Given that, as our main result we have esti- 
mated the age of the family. Using two different methods, MCC 
(Tsiganis et al. 2007; Novakovic et al. 2010a) and backward inte- 
grations ( Nesvornv et al.ll2002h . we estimate the age of the family 
to be 6.9 ± 2.3 Myr. Thus, we establish it as another young asteroid 
family. 

The paper is organized as follows: In Section |2] we present 
the basic properties of the Theobalda family. We use the hierarchi - 
cal clustering method (HCM), proposed by IZappala et al.l ( Il99(]r) . 
to identify family members. Next, the cumulative absolute mag- 
nitude distribution of identified family members is discussed, and 
the size of parent body is estimated. The dynamical characteris- 
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tics of the region occupied by Theobalda family are presented and 
discussed, and main mean motion and secular resonances, in that 
region, identified. The dynamical stability of the family members 
is analyzed, in particular the stability of the largest member of the 
family (778) Theobalda. In Section [3] we estimate the age of the 
family. This is performed firstly by using the backward integration 
method, and then by using the method of chaotic chronology. The 
good agreement between these two results indicates a reliable age 
determination. Finally, in Section|4]we summarize our results, dis- 
cus some possibly interesting relations to other works, and draw 
our conclusions. 



2 THEOBALDA FAMILY: THE BASIC FACTS 

This asteroid family has attracted little attention so far, mainly 
because the number of asteroids associated with it was relatively 
small. However, the situation is different at present, and, as we will 
show later, Theobalda family now has over 100 known members. 
This number is large enough that the family characteristics can be 
reliably determined. 



2.1 Identification of tlie family members 

The identification of family members is the first step in our study 
of the Theobalda family. This is done by applying the HCM to 
the catalog of synthetic proper elements of numbered asteroids 
JKnezevic & M ilani 2000, l2003h from AstDyf] (database as of Oc- 
tober 2009). The HCM requires that distances among the family 
members, in the proper elements space, are less than the so called 
cut-off distance (dcutoff), which has dimension of velocity. As 
the cut-off distance is a free parameter of HCM, we tested differ- 
ent values ranging from 20 to 135 ms^^. Also, we apply HCM 
using two different central objects: (i) (778) Theobalda which 
has a chaotic orbit, and (ii) (84892) 2003(5-D79 which is on the 
relatively stable orbit. The results are shown in the top panel of 
Fig. [T] The HCM identified the family around (778) Theobalda 
for dcuto// > 60ms-\ while around (84892) 2003QI>79 fam- 
ily exists even for lowest tested value of dcutoff ~ 20 ms~^ . For 
dcutoff > 60 ms^^ resulting family is the same. This suggests 
that (778) Theobalda probably has been displaced from its origi- 
nal position due to the chaotic diffusion[j In the bottom panel of 
Fig. [T]the best-fit power-law index 7 of the form N(<H) a 10^^ 
of the cumulative absolute magnitude (H) distribution in the range 
HG[I3-15], as a function of cut-off distance (dcutoff), is showrQ 
For dcutoff €[75,115] ms~^ the number of asteroids as well as in- 
dex 7 are nearly constant, and probably each value from this inter- 
val can be safely used to identified family members by HCM|j We 
adopted value of dcutoff ~ 85 ms^^ to identify nominal family. 
For this value of dcutoff, HCM linked 128 asteroids to Theobalda 
family. There are two main reasons for our choice. The first one 



^ http://hamilton.dm.unipi.it/astdys/ 

^ Note, that this is very similar to Veritas family, a nd situation with the 
largest member of this family (490) Veritas jTsiganis et al. 2007). 
^ Instead of the index 7, the exponent of the cumulative distribution 
can be obtained in terms of diameters rather than absolute magnitudes 
jDeirOro & CellindllOOTl) . However, as we do not know albedos for most 
of the asteroids, necessary to convert from absolute magnitudes to diame- 
ters, we chose to work with 7. 

* Usually one adopts the value of dcutoff that corresponds to the center of 
the interval over which the index 7 is constant jVokrouhlickv et alj|200q) . 
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Figure 1. Number of asteroids associated with Theobalda family (top) 
and a power-law index 7 of the cumulative magnitude distribution in the 
range HG[13-15] (bottom) as a function of the cut-off distance (dcutoff)- 
Note that in both cases the respective values are nearly constant for 
dcutoff G[75,95] ms~^. In the bottom panel horizontal solid line and 
two dashed lines represent background level and its eiTor-bars respectively. 
In the top panel two obvious critical values of dcutoff ^re: (i) 60 ms~^, 
when group around (778) Theobalda merge with the group around (84892) 
2003QDrQ, and (ii) 130 ms~^, when the family starts to merge with the 
local background population. However, the number of asteroids associated 
with the family is constant for dcutoff S [75,95] ms~^. For the nominal 
family we chose dcutoff = 85 ms~^ (see text for additional explanation). 



is that this value of velocity cutoff corresponds to the center of the 
plateau which can be seen in Fig.[T] The second reason is very good 
agreement between the ages of family estimated applying MCC to 
two different groups of family members. We will explain this in 
more detail in Section[3] 

Note that the values of 7 for family members are always 
larger than the value of the same index calculated for background 
asteroids. This is the first indication that the family is relatively 
young. On the other hand, Parker et al. (2008) estimated 7=0.44 for 
HG[13.0-15.5]. This value is much lower than ours, as we found 
7=0.60±0.02 for the nominal family, and very close to the value 
that we found for background popula tion in the region of Theobalda 
family. Probablv. lParker et al.l ( l2008h underestimated this value due 
to the observational incompleteness, as they worked in the range 
HG[13-15.5] and used smaller dataset for which Sloan Digital Sky 
Survey (SDSS) colors were available. Although they linked 100 as- 
teroids with the family, a significant number of these asteroids are 
probably interlopers. 
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Figure 2. DisMbution (top) and cumulative distribution (bottom) of the 
members of Theobalda asteroid family as a function of the absolute magni- 
tude. 



2.2 Size of the parent body 

To estimate the diameter of the parent body (Dpb) of the 
Theobalda family, it is necessary to account for small and still 
undiscovered family members. In gen eral, data set on astero ids be- 
low H=15 mag is basically complete JGladman et al.ll2009l) . How- 
ever, as we are dealing with the family at the edge of outer belt, the 
family members are the C-type asteroids which are several times 
darker than the S-type asteroids, and, since we are using catalog of 
synthetic proper elements which does not include all known aster- 
oids, the completeness limit for our sample has to be analyzed. An 
indication about the completeness limit can be obtained by simply 
looking at Fig. |2] for the value of H where slopes of two distribu- 
tion curves change. This is approximately at about 7^=14.5 mag. 
Somewhat better estimation of completeness limit can be inferred 
using the catalog of asteroids, which are not included in the catalog 
of synthetic proper elements we deal with, i.e. the catalog of multi- 
opposition objects maintained at AstDys web site. As about 99 per 
cent of multi-opposition objects, with osculating semi-major axes 
in the range [3.15,3.20] au, have H >14.2 mag (see Fig. [3), we 
assume that the catalog of synthetic proper elements of asteroids in 
this region is complete up to H=14.2 mag. 

In order to overcome the problem of observational incom- 
pleteness, using the size-frequency distribution (SFD) of the main 
belt asteroids estimated by Gladman et al. ( 2009 jj and the fact that 
SFD of asteroid families is consid ered to be som ewhat shallower 
than that of the background CMorbidelli et al.l2003f) we added some 
fictitious bodies with H £[14.2,17.0] to the family. This was per- 
formed in such a way to make SFD of "extended" family (real aster- 
oids -I- fictitious objects), to be somewhat shallower than obtained 
SFD for background asteroids. More precisely, in order to be able 
to estimate the uncertainty of our approach, we generated 100 dif- 
ferent sets of fictitious objects and then estimated the size of the 
parent body from each of the 100 sets. 

DpB corresponds to a spherical body with volume equal 



^ The asteroid size distribution at diameters D<10 km is still poorly 
known. Various models and extrapolations yield very different esti- 
mates of the number of km-sized and smaller main belt asteroids. How- 
ever, other possible estimates (e.g. Ivezic et al.] |200ll : lTedesco et alJIlOOSi ; 
Iwiegert et al. 2007; Yo shida & Nakamura .20071) would not affect our re- 
sults significantly. 
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Figure 3. The distribution of multi-opposition asteroids, in the (a,H) plane, 
from AstDys web site (database as of October 2009). The catalog consists 
mostly of the objects discovered more recently than objects included in the 
catalog of synthetic proper elements. Thus, it provides good opportunity 
to estimate completeness limit of the catalog of synthetic proper elements, 
which is marked by a red line. The green dots represent the members of 
Theobalda family. 



to the estimated total volume of all the family members, includ- 
ing these added members with H G[14.2,17.0]. Next, we as- 
sume all family members have the same geometri c albedo (p„) 
as (77 8) Theobalda, that is p^=0.0589 according to iTedesco et al.l 
1 120021) . Having the values of H and Vy, the radius i? of a body can 
be estimated, using the relation (e.g. lBowell et aljl989i) 



R (km) = 1329 



10" 



/Pv 



(1) 



This allows us to infer that the diameteu of the parent body was 
DpB « 78 ± 9 km. The estimated uncertainty accounts for uncer- 
tainties in albedos, absolute magnitudes and SFD. Also it takes into 
account the dependence on the HCM cut-off. However, the real un- 
certainty is somewhat larger, e.g. because of the possible interlop- 
ers. According to our estimate, the largest remnant (778) Theobalda 
contains 87 per cent of the mass of the parent body. Although this 
should be considered as an upper bound, because it does not ac- 
count for small family members with H > 17, our result suggest 
that Theobalda family was produced by a cratering impact. The 
typical density ( p) for C-type asteroids is p = 1500 kg m"'' (e.g. 
iBroz et al.ll2005f) . Given that, the escape velocit}[j from Theobalda 
family parent body was Vesc~ 32 ms~^. 



2.3 Dynamical characteristics 

The dynamics in the region of the phase space occupied by 
Theobalda family members is much like in the case of the Veri- 
tas family because two families stretch over the similar range of 



^ Probably the better way to es timate the size of the parent body is the one 
proposed bv lDurda et alJJ2007D . However, we were unable to find appropri- 
ate match with SFDs pubUshed in that paper using simple visual comparison 
of plots. 

^ Compensating for collective effects in the cloud of dispersing fragments, 
Vesc = 1.64 X GM/R, where GM is the product of the gravitational con- 
stant and the parent b ody mass, R is radius of the parent body, while 1 .64 
is an empirical factor IPetit & Farinellall 19931) . 
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Figure 4. The dynamical structure of the region occupied by Theobalda 
family along with surrounding area. The color scale codes Lyapounov Char- 
acteristic Exponents (in the units of 10^® yf'^^) for 10,000 test pailicles. 
The ellipse represents assumed positions of the Theobalda family members 
immediately after break-up. 



the proper semi-major axes. The dynam ics in the region of Veri- 
tas family is very well stud i ed (see e.g. Milani & Farinella 1994; 
Knezevic & Pavlovid l2002l : iNesvomv et al.. .2003. : .Tsiganis et al.. 



20070 . Therefore, here we will focus only on some differences be- 



tween dynamics in the regions occupied by two families. The dif- 
ferences arise from the fact that proper eccentricities of Theobalda 
family members (ep~0.25) are significantly higher than those of 
Veritas family (ep~0.06). Also, the proper inclinations are by 
about 5° higher. These make Theobalda family members even more 
strongly chaotic than Veritas family members. 

In Fig.|4]the Lyapounov Characteristic Exponents (LCEs), as 
a function of proper semi-major axis and eccentricity, in the region 
occupied by Theobalda family members, as well as in the surround- 
ing area are shown. Most of the orbits in that region are unstable 
even for comparatively low values of proper eccentricity, while for 
eccentricity above 0.3 almost all chaotic zones are connected form- 
ing a wide chaotic sea with a fast diffusion therein. The chaos is 
also dominating in the region where Theobalda family is located 
(inside or close to the equivelocity ellipse). This can be better ap- 
preciated from Fig. [5] where LCEs of Theobalda family members 
are shown. The vertical strip of the largest at ag~3. 174 au values of 
LCEs is associated with (5, -2, -2) three-bod5lj MMR, but it seems 
that this chaotic zone includes (3, 3, -2) and (7, -7, -2) three-body 
MMRs as well. Most of the bodies have LCE >lxlO~'' yr'^, 
which corresponds to the Lyapounov times Tiyap < 10,000 yr 
and these bodies are probably in t he so-called Chirikov regime 
JGuzzo et alj2002l : lMorbidelll2002l) . 

Using the proper frequencies g (average rate of the perihelion 
longitude -go) and s (average rate of the node longitude f2Jj we 
found that Theobalda family region is also crossed by two secular 
resonances g + s — (75 — se and <? + s — ge — se (Fig.|6j. Both are 
of the order 4, i.e. they arise from the pe rturbing terms of degree of 
at le ast 4 in eccentricity and inclination (Milan i & Knezeviall990l : 
iKnezevic et al.il991i) . However, we did not find evidence (see Sec- 
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Figure 5. The same as in Fig.|4]but for LCEs of the real family members. 
The linear interpolation is used in order to cover the complete region shown 
in the figure. 



tion l3.2.l1 l that these resonances increase diffusion speed. Probably 
this is because these resonances are effective only in narrow bands 
within the Theobalda familv^°l Consequently, some of the family 
members, during their secular cycles, might be temporally trapped 
in one or both secular resonances, but most of the time these as- 
teroids are outside the secular resonances. In Fig. [6] the time evo- 
lution of the critical angles (Ti=cc7-l-f2-cc75-f26, and a2='cu+Q,-zuQ- 
Q.(i, for asteroid (778) Theobalda are shown. This asteroid might be 
temporally trapped in both secular resonances. Although the short 
episodes of "libration" are visible, these events may be related to 
resonance crossing rather then to the resonance trapping. Most of 
the time both critical arguments circulate. 

In Fig. [7] distributions of family members, as identified by 
HCM, are shown along with the positions of main mean motion 
and secular resonances. Obviously, the structure of the family is a 
result of dynamical mechanisms at work, which are mainly con- 
trolled by MMRs. The largest spread of family members, in both 
(ap,ep) and (ap,Ip) planes, is associated to (5, -2, -2) resonance. 
Somewhat smaller spread is observable in the (3, 3, -2) resonance, 
while (7, -7, -2) resonance caused only small diffusion of asteroids. 
This agree very well with obtained values of LCEs (see Fig. [5). 

It is interesting to note that there are gaps (without family 
members) between the (3, 3, -2) and (5, -2, -2) resonances, as well 
as between the (5, -2, -2) and (7, -7, -2) resonances. We suggest 
that this is another confirmation that all these three resonances are 
connected and make one wide chaotic zone. Because of this, all as- 
teroids from ap~3.167 au to ap~3.181 au reside in one of these 
three resonances. The asteroids can switch from one resonances 
to another (see Fig. |9j, but on the time scale of a few Myr this 
is a rare event, so that each asteroid spends most of the time in 
one of the resonances. As a result, due to the some uncertainty 
in the procedure of computation of synthetic proper elements for 
resonant asteroids, i.e. averaging does not work well, all bodies 
appear to be located in (or close to) the center of one of the res- 
onances. T his can be verified by usin g the analytical proper ele- 
ment^ of iMilani & Knezevig ( Il99(]|) . These elements are calcu- 



* All three-body mean motion resonances jNesvomv & Morbidellilll998l) 
discussed in this paper are among Jupiter, Saturn and asteroid. 
^ The secular freq uencies of Jupiter (g^) and Saturn (36. ■56) are taken from 
lNobilietalJil989 



The secular resonance g - 



-56- 



sg has much more influence on dy- 






namics of family members in the case of Padua family (see Carruba 200' 
^^ We did not use analytical proper elements, in our other analysis through- 
out the paper, because they are not enough accurate in this high eccentricity 
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Figure 6. Time evolution of the critical angles o-i=iAj+Q--n:'5-Qc, (top) 
and a-2=iAj+Q-zu(;-Qfi (bottom) of the secular resonances for period from 
1 .2 Myr to 1 .6 Myr in the past. The short episodes of trapping inside the res- 
onance seem possible, but several events of resonance crossing are clearly 
visible. Reversal of direction of circulation is related to periods when the 
orbit interact with these secular resonances. 



lated by means of analytical theory based on the series develop- 
ment of the perturbing Hamiltonian, and which does not include 
averaging. The distribution of Theobalda family members in the 
space of analytical proper elements (in the ap, Cp plane) is shown 
in Fig. [S] The shown distribution is roughly random and without 
gaps in terms of proper semi-major axis, what confirms our claim 
that the gaps appear due to the averaging procedure. Moreover, it 
means that switching from one resonance to another must be a rare 
event, but, the fact that not all of the asteroids are located in the 
center of one of the resonances, is another evidence that resonance 
switching is possible, i.e. these three resonances are connected. 

The position of the largest remnant, asteroid (778) Theobalda, 
is not close to the center of the family. This is evident also in 
(Upjp) plane, but it is more obvious in (ap,ep) plane. As we al- 
ready mentioned above, this asteroid has probably been displaced 
from its original position due to the some dynamical mechanisms. 
It is located close to or inside the (7, -7, -2) three-body MMR, 
which might be responsible for its relatively high proper eccen- 
tricity. However, its proper semi-major axis is also larger than that 
of the center of family, and this could not be explained by (7,-7, -2) 
resonance. Because of that, we investigate dynamics of this aster- 
oid in more detail. The orbit of (778) Theobalda is propagated for 
100 Myr back in timg^H As we will show later, the family is about 
6-7 Myr old. Why than it is meaningful to integrate 100 Myr? The 



region. This can be appreciated comparing the distributions of Theobalda 
family members with ap>3.183 au, shown for two different kinds of proper 
elements. Obvious grouping of the regular members, which is clearly vis- 
ible in the space of synthetic proper elements (Fig. |7), disappears in the 
space of analytical proper elements (Fig. [8). 

^■^ All integrations presented in this paper are performed using the public 
domain ORBIT9 integrator embedded in the multipurpose OrbFit package 
(http://hamilton.dm.unipi.it/astdys/ 1, and dynamical model that includes the 
four major planets (from Jupiter to Neptune) as perturbing bodies. The in- 
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Figure 7. Distribution of the known Theobalda family members in the 
(ap,Ip) plane (top) and (ap,ep) plane (bottom). The superimposed el- 
lipses represen t equivelocity curves, c omputed according to the equations 
of Gauss (e.g. iMorbidelli et al.lll995l) . for a velocities of ii = 35 ms~^ 
(inner) and d = 40 ms~^ (outer), true anomaly / = 85° and argument 
of pericentre uj = 95°. The blue points represent family members identi- 
fied for d^^taf f = 65 ms~^, while red points represent additional fam- 
ily members identified for d^^i^jj = 85 ms~^. The size of each point 
corresponds to the diameter of the body. The green dashed lines mark ap- 
proximately borders of three-body MMRs, while pink dashed lines show 
locations of secular resonances. 
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Figure 8. The same as in Fig.|2](bottom), but for analytical proper elements. 
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Figure 9. Time evolution for 100 Myr back-in-time of the mean semi-major 
axis for asteroid (778) Theobalda. For the first about 15 Myr of integration it 
resides in (7, -7, -2) resonance, but then, it switches to (5, -2, -2) resonance 
where it spends most of the time covered by integrations. Also, it exhibits 
short episode of trapping inside (3, 3, -2) resonance at about 42 Myr. 



answer is hidden in the chaotic motion of this asteroid. As we know, 
chaos is not predictable on the time scales of several times the in- 
verse of LCE, which is in the case of (778) Theobalda ~8,000 yr. 
All that we can achieve is to show what kind of behavior (i.e. mo- 
tion) is possible. In this respect, our 100 Myr long integrations are 



equivalent to many shorter integrations wit h slightly differ ent initial 
conditions. Similar technique was used bv lLaskan ( Il994h to study 
stability of the Solar system. 

Fig. |9] shows 100 Myr of the back-in-time evolution of the 
mean semi-major axis of asteroid (778) Theobalda. Initially, its 
semi-major axis oscillates around 3.181 au (close to (7, -7, -2) 
resonance), but after about 15 Myr (in the past) the value of the 
semi-major axis drops to ~3.174 au. Around this time it actually 
switches from (7, -7, -2) to (5, -2, -2) resonance. There is another 
switch of the resonance about 42 Myr, when the asteroid is tem- 
porally trapped in (3, 3, -2) resonance at ~3.168 au. Finally, at 
about 95 Myr in the past, it went back from (5, -2, -2) to (7, -7, -2) 
resonance. This is the confirmation that these three resonances are 
connected at higher eccentricities (e^ > 0.25). But again, due to 
the chaoticity, the behavior of the mean semi-major axis of aster- 
oid (778) Theobalda does not represent quantitatively its real mo- 
tion, but qualitatively. Still, behavior of its semi-major axis suggests 
that chaos may be responsible for displacement of (778) Theobalda 
from the center of family in terms of ep and Ip. However, if this 
was the case, this asteroid probably spend some time residing in 
(5, -2, -2) resonance which is strong enough to increase its eccen- 
tricity from ep«0.253 to ep«0.259, on the time scale of several 
Myr. 

Studying the distribution of family members shown in Fig.|7] 
it can be noted that there are no family members located inside 
equivelocity ellipses, at ap~3.165 au. Contrary to the gaps between 
the resonances, the absence of asteroids that belong to the family 
in this region cannot be explained by dynamical instability or by 
"weakness" of the procedure of proper elements calculation. Al- 
though, a detail study of this problem is beyond the scope of our 
work, we believe that this may be related to the impact character- 



direct effect of the inner planets is accounted for by applying a barycentric 
correction to the initial conditions. 



istics (cratering event), which "forced" fragments to be symmetri- 
cally distributed around the semi-major axis of the largest fragment 
(778) Theobalda. 



3 THE AGE OF THEOBALDA FAMILY 

3.1 Backward integration 

Backward integration of orbits is very accurate method for family 
age estimation, which works well with young families. It is based 
on the fact that due to the planetary perturbations the orientation 
of orbits in the space changes over time. Consequently, two angles 
that determine the orientation of orbits in space, the longitude of 
the ascending node (Q) and the longitude of perihelion (zu), evolve 
with different but nearly constant speeds for individual orbits. Af- 
ter some time this effect spreads out Q and zu uniformly around 
360°. On the other hand, immediately after the disruption of the 
parent body, the orientations of the fragments' orbits must have 
been nearly the same. Given that, the age of an asteroid family 
can be determined by integrating the orbits of the family members 
backwards, until the orbital orient ation angles cluster around so me 
value. This method was used by iNesvornv et al.l ( 120021 l2003h to 
determine the ages of the Karin cluster (5.8±0.2 Myr) and Veritas 
family (8.3±0.5 Myr). 

Here we applied Nesvomy et al.'s method to try to estimate the 
age of Theobalda family. By integrating the orbits of the Theobalda 
asteroid family back in time, hopefully, we can find a conjunction 
of orbital elements (fi and zu), which occurred only immediately 
after the disruption of the parent body. This method is, however, 
limited to groups of objects moving on regular orbits, which, even 
in that case, can be accurately track up t o 20 Myr in the past . 
Similarly as in the case of Veritas family JNesvomv et al.ll2003l : 
iTsiganis et al.ll2007r) . only a fraction of Theobalda family mem- 
bers satisfies this condition and can b e accurately integrated back 
in ti me. Also, as was poi nted out bv INesvornv et al.l J2003h (see 
also Nesvomy et al.ll2008r) the region around ap=3.175 au is close 
enough to the 2/1 MMR with Jupiter to undergo fast differential 
evolution of the arguments of perihelion. This induces variabil- 
ity in the evolution histories and complicates any attempt to de- 
termine the age of the Theobalda family using arguments of per- 
ihelion. Thus, we selected 13 Theobalda family members which 
have Lyapounov times Tiyap > lO^'yr and propagated their orbits 
20 Myr backwards. All these members are located at ap>3. 183 au. 
In Fig. [To] the average value of Afl, for these 13 asteroids, is 
shown. Conjunction of nodal longitudes at ~6.2 Myr suggests that 
the Theobalda family, or at least a part of the family located at 
ap>3.183 au, was formed by a catastrophic collision at that time. 
The average Af2, at ~6.2 Myr, is ~58°, much smaller than at any 
other time. This suggests a statistical significance of the ~6.2 Myr 
event. In this case, however, (ASl) values are substantially more 
spread at ~6.2 Myr than in the case of Karin cluster ((AS!) is 
~10°) or Veritas family (( ASl) ~40°). This is primarily due to two 
reasons: (1) at least a few MMRs exist in the semi-major axis range 
from 3.18 to 3.19 au; thus, despite the present long Lyapounov 
times of the selected orbits, these orbits might have experienced 
periods of chaotic motion in the past; and (2) all regular bodies, 
whose orbits can be accurately tracked back in time, are small bod- 
ies (<5 km) and consequently subject to Yarkovsky thermal force, 
which, even on this relatively short time scale, can produce large 
enough changes in the semi-major axes, and consequently to affect 
the secular frequencies in a way that is difficult to reconstruct. 
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In order to estimate how sensitive this result is on the semi- 
major axis drift due to the Yarkovsky effect, an additional investiga- 
tion should be carried out. As the Yarkovsky induced drift depends 
on several parameters, we had to decide the values of the parame- 
ters characterizing it. These are asteroid spin axis orientation (7), 
rotational period (P), surface and bulk densities (p), surface ther- 
mal conductivity (K) and specific heat capacity (C). As Theobalda 
family members are most likely C-type asteroids, we have adopted 
thefollowing values of parameters: K = 0.01 — 0.5 [W(m K)^^], 
C = 680 — 1500 [J (Kkg)~^], and the same value for surface and 
bulk density p = 1300 — 1500 [kgm~^]. The rotational periods 
are chosen according to a Gaussian distribution peaked at P = 8 h, 
while the distribution of spin axes orientation is assume d to be uni- 
form. These values are cons istent with C-type asteroids ( lBroa2006l : 
iBroz & Vokrouhlickvll2008l) . 

Next, we made 20 "yarko" clones for each of the 13 regu- 
lar members, by assigning random values of the parameters, from 
adopted intervals, to each clone. Then, we integrateco the orbits 
of all clones (260 orbits in total), but accounting not only for grav- 
itational perturbations, but also for Yarkovsky effect. The initial 
orbital elements of the asteroids and planets were the same as in 
the previous experiment. Finally, we checked how the value of av- 
erage ASl change with different combinations of clones. We found 
that the result shown in Fig.[TO]is very sensitive to the Yarkovsky 
induced drift, as expected. In a few cases any significant clustering 
even disappeared, but in most of the cases we obtained a deeper 
minimum. The deepest minimum that we found is related to the 
clustering within «31° at about 6.4 Myr ago (Fig. llOt , which 
is still within the error bars obtained from the integrations with- 
out Yarkovsky force. We would like to note here that this high 
sensitivity of the result on the Yarkovsky parameters could help 
us to estimate the rotational periods and spin axis orientations of 
these 13 asteroids. This ca n be achieved similarly as was done by 
iNesvornv & Bottkd ( 12004) for the Karin cluster members, but we 
reserve this for a future work. 

Although, the clustering at about 6.2 Myr within ~58° is the 
most significant on the time scale of 20 Myr, there is another clus- 
tering at about 15.5 Myr within ~65° (see Fig. llOt. As this clus- 
tering appears in the more distant past where we should expect less 
tight clustering, it is not possible to rule out its significance. Also, to 
use the argument of perihelions is impossible, because the changes 
in the semi-major axes, caused by Yarkovsky effect, coupled with 
large gradienO of secular frequency (dg/da ~0.3° yr^^ au, 
where g is the longitude of perihelion frequency), erase evolution 
histories of these angles. Given that, we believe that, in the case of 
the Theobalda family, backward integration method is not enough 
to draw a firm conclusion about the age of the family. 



3.2 Chaotic chronology 

In this section we present results obtained by using MCC in order to 
estimate t he age of Theobalda fami ly. This model was successfully 
applied bv lNovakovic et al.l ( l2010ah to estimate ages of Veritas and 
Lixiaohua asteroid families (see also iNovakovic et al. 2010b) . In 
order to apply MCC the family has to be located in the region of the 




Figure 10. The average differences in nodal longitudes ((An)) for 13 
members of Theobalda family, with regular orbits. The results obtained in a 
purely gravitational model (top) and with yarko clones (bottom), are shown. 
The most important feature is clustering at about 6.2 Myr ago within f»58°. 
This minimum becomes significantly deeper (f»31°) with yarko clones, and 
also slightly shifted to 6.4 Myr ago. 



main belt where diffusion takes place. Also, it is necessary that dif- 
fusion is fast enough to cause measurable effects, but slow enough 
so that most of the family members are still forming a robust family 
structure. As our results about diffusion speed suggest, Theobalda 
family is an excellent example in this respect (see Section [3.2.1b . 

The basic steps and the model which we used in our 
Markov Chain Monte C arlo (MCMC) simulations are explained in 
INovakovic et al.l feOlOah . and thus we will describe these here only 
briefly. Our model simulate the evolution of the family in the 3-D 
space, i.e. proper semi-major axis Op and two actions Ji, J2 (see 
Section [3.2.1l for definition of these actions). At the beginning of a 
simulation the random walkers are distributed in the region which 
was presumably occupied by the family members immediately af- 
ter the impact event. Then, at each time step dt the random walkers 
can change their positions in every direction, in the 3-D space. The 
length of the jump in Sp is contr olled by Yarkovsky thermal force 
JFarinella & Vokrouhlickvll 19991) . while the length of the jumps in 
Ji and J2 depend on diffusion speed, i.e. on the diffusion coeffi- 
cients. At the time step when 0.3 per cent of random walkers leave 
an ellipsq^in the (Ji, J2) plane, which corresponds to a 3a confi- 
dence interval of a two-dimensional Gaussian distribution, the sim- 
ulation stops. The number of time steps multiplied by the time step 
dt gives the age of the family. 



^^ These int egrations were performe d using 0RBIT9 integrator in the Grid 
environment INovakovic et alj2009l) . 

^* Caused by the proximity of the Theobalda family to the 2/1 resonance 
with Jupiter 



^^ The ellipse is determined by the present size of the family or, as in this 
case, by the present size of particular part of the family. It should not be 
confused with equivelocity ellipses shown e.g. in Fig.]?] 
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be noted here that this might affect age estimation using backward 
integration method, but not significantly, because the most of 13 as- 
teroids that we used to apply backward integration method are lo- 
cated close to a^ =3.185 au, where both values, D{Ji) andD(J2), 
are close to zero. 

The important general conclusion can be drawn by comparing 
the values of diffus ion coefficients obtained for the region occupied 
by Veritas family JNovakovic et alJuOlOah to the values obtained 
here. Two families stretch over the similar range of the semi-major 
axis, but members of Theobalda family have somewhat higher in- 
clinations and significantly higher eccentricities. The estimated dif- 
fusion is about one order of magnitude faster in the region occupied 
by Theobalda family than in the region occupied by Veritas family. 
This confirms the fact that chaos is dominant at higher eccentrici- 
ties. 



Figure 11. The values of diffusion coefficients D{J\) (red) and D{J2) 
(blue) in the Theobalda family region, shown here as functions of the proper 
semi-major axis ap. Note that D{J2) is practically zero for ap >3.18 au. 



3.2.1 Dijfusion coefficients 

One of the most important information, which are needed as input 
for MCMC simulations, are the values of diffusion coefficients i n 
the region of interest. As was shown bv lNovakovic et al.l j2010ah . 
to obtain good estimate of the family age by MCC, it is enough to 
determine diffusion coefficients as a function of proper semi-major 
axis ap. This is our next step. 

As well as MCC, the procedure of determination of diffu- 
sion coeffici ents, as the func t ions of proper semi-major axis, is 
described in iNovakovic et al.l ( 1201 Oaf) . Let us mention here only 
its main features and numbers related to this work: the orbits of 
~5,000 fictitious bodies distributed randomly in the same ranges of 
osculating orbital elements as the real family members at present, 
are pr opagated for 10 Myr; th en, the time series of mean ele- 
ments JMilani & Knezeviall998 ) for all of them are calculated; the 
mean elements are transformed to actions according to relationa^^l 



Ji 



and J2 



■ sin Im. ; next, the family is 



split in the small cells, in terms of ap, using a kind of moving- 
average technique with cell size of Aop = 5 x 10"'' au and step 
size of Sap = 2 x 10~ au; finally, the mean squared displacements 
{(AJ)^), for both actions, are calculated, and the diffusion coeffi- 
cients Z)(Ji), D(J2) for each cell as the least-squares fit slope of 
the {(A J) ) (t) curve, are determined. 

The obtained values of diffusion coefficients -D(Ji) and 
D{J2), in the Theobalda family region are shown in Fig.[TT] The 
fastest diffusion is associated to (5, -2, -2) three-body MMR, but 
the diffusion is very fast in (3, 3, -2) three-body MMR as well, and 
these two chaotic zones seem to be connected. The third chaotic 
zone, associated to (7, -7, -2) resonance, is connected to the first 
two in terms of diffusion in Ji, but not in terms of diffusion in 
J2. This is in a relatively good agreement with results presented 
in Section [23] The diffusion is somewhat faster in Ji than in J2, 
and while the local minimum near the center of (5, -2, -2) reso- 
nance exist fo r D(J\), simi larly as in the case of Veritas family 
(see Novakovic et al.ll2010ah . while there is no such a feature for 
D{J2). In the region for ap >3.18 au the values of D{J2) are 
practically zero, but there is some diffusion in the Ji. It should 



° In these relations aj denotes Jupiter's semi-major axis, em the mean 
eccentricity and Im the mean inclination of the asteroids. 



3.2.2 Monte Carlo simulations 

Having obtained the values of diffusion coefficients, we are ready 
to apply MCC to estimate the age of the family. There are two sepa- 
rate parts of the Theobalda family suitable for application of MCC. 
These are bodies inside (5, -2, -2 ) and (3, 3, -2) three-body MMRs. 
Following iTsiganis et alj J2007D who deal with Veritas family, we 
called these bodies Group A (5, -2, -2) and Group B (3, 3, -2). As 
the results about diffusion coefficients confirmed, there exists sig- 
nificant diffusion in both groups. This gives an unique opportunity 
to apply MCC to both groups and to obtain two independent age 
estimates. A good agreement between these two estimates, as well 
as with the age derived using backward integration method, would 
suggest a reliable result. 

As the present size of the chaotic zone is a critical parameter in 
our model, we start with the family as identified by applying HCM 
for velocity cutoff of dcutoff ~ 65 ms~^. This is probably the 
lowest acceptable value of dcutoff in the case of Theobalda family. 
With this cutoff velocity we identified 30 bodies from Group A and 
16 bodies from Group B. Corresponding sizes of these groups in Ji 
and J2 are: Group A (10.67 ± 1.03) x 10"'' and (3.82 ± 0.47) x 
10"^ Group B (10.32 ± 1.05) x 10"* and (4.00 ± 1.06) x 10""*. 

Using these sizes of two chaotic groups, for each group, we 
performed 16 sets of MCMC simulations (each set consisting of 
100 runs), by using different number of random walkers n (2000 
or 5000), time step dt (from 100 yr to 2000 yr) and for two initial 
sizes of the family which correspond to velocities of u = 35 ms"'^ 
and n = 40 ms"^ (see Fig.[7]l. From these simulations we derived 
the age of family to be 2.5 ± 1.1 Myr (using Group A bodies) 
and 7.2 ±3.1 Myr (using Group B bodies)|^ The obvious dis- 
crepancy between two results needs to be investigated further. The 
age obtained from Group B is in agreement with what we found 
using backward integration method, while the age obtained from 
Group A suggests that the family could be much younger. Also, the 
discrepancy between ages derived from two different groups, may 
be an indication that identification of family members has not been 
good, i.e. the velocity cut-off of 65 ms^'^ is too low. 

Because of that, we repeated all simulations using our nom- 
inal family. For dcutoff ~ 85 ms"^ we identified 40 bodies 
from Group A and 18 bodies from Group B. In this case, the 
corresponding sizes of these groups in Ji and J2 are: Group A 



^^ The main source of the eiTor is uncertainty in the determination of the 
present size of the Group B, due to the small number of members in this 
Group. 
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Figure 12. The age of Theobalda family derived using different methods, 
asteroids and parameters: (i) MCC applied to Group A for d^^^^fj-- = 
85 ms~^; (ii) MCC applied to Group A for cicuto// = 65 ms~^; (iii) 
MCC applied to Group B for dcutoff = 85 ms~^; (iv) MCC applied to 
Group B for d^^t^ff = 65 ms~^; (v) age obtained by backward integra- 
tion method. The bold horizontal line and dashed area correspond to our 
final estimate of Theobalda family age and its eiTor respectively. 



(14.86 ± 1.24) X 10"" and (6.03 ± 0.65) x 10"*; Group B 
(11.87 ± 1.29) X 10"* and (3.80 ± 0.95) x 10"*. From these 
sizes and the same sets of MCMC simulations as in the previous 
case, we derived the age of Theobalda family to be 6.9 ±1.8 Myr 
(using Group A bodies) and 7.2 ± 3.0 Myr (using Group B bodies). 

Now, the agreement between two results is very good, and 
also, both results agree quite well with the age obtained by back- 
ward integration methodic This, in our opinion, is a very strong 
indication that Theobalda family was formed about 7 Myr ago. 

The fact, that four out of five, different age estimates, agree 
well, is the reason why we reject the age of 2.5 ± 1.1 Myr, de- 
rived using Group A (dcutoff = 65 ms~^), as a possible solution. 
Thus, in order to obtain our final estimate of the age of Theobalda 
family, we use four results which are in a good agreement (see 
Fig.ll2t. The values of the mean (/x) and standard deviation (a) of 
non-overlapping sub-samples, of the same size, can be calculated 
as: 

/. = ^^^i^ (2) 



^ZMr- 



l)crf + rfij) — rnr^j? 



(3) 



where ra is the number of sub-samples, r is the size of each sub- 
sample (100 in our case), /i; is the mean of i-th sub-sample, and 
(T,; is the standard deviation of i-th sub-sample. Using Eqs.|2]and 
[3] we obtain the final age estimate, Theobalda asteroid family is 
6.9 ± 2.3 Myr old. 



^* The good agreement between ages obtained applying MCC to groups A 
i B, is one of the reasons why we adopted value of di-^t^ff = 85 ms^^ 
to identify nominal family. It should be noted here that identification of 
resonant family members is not straightforward. Too small cut-off, on one 
hand, may prevents identification some of the real family members. On the 
other hand, the large cut-off could associate some interlopers with family. 
This may be one of the reasons for the variations of age estimates obtained 
by MCC for two different cut-off velocities. 



4 SUMMARY, DISCUSSIONS AND CONCLUSIONS 

We have presented here a detailed study of Theobalda asteroid fam- 
ily. We found that family now consists of 128 members. By ana- 
lyzing SFD of the identified family members we were able to in- 
fer diameter of the parent body to be Dps « 78 ± 9 km. How- 
ever, this estimate is based on the assumption that all family mem- 
bers have the same albedo as the largest family member, aster- 
oid (778) Theobalda. In order to obtain better estimate, the albe- 
dos of as many as possible family members are desirable. Ongo- 
ing projects, such as Wide-Field Infrared Survey Explorer (WISE) 
should improve situation significantly in this respect. 

The most, but not all, of Theobalda family members move on 
chaotic orbits, thus, giving rise to the significant chaotic diffusion 
which has been changing the kinematical structure of the family 
over time. The study of dynamical characteristics, in the region oc- 
cupied by the family, showed that three three-body MMRs are the 
most efficient in shaping the family. These are (3, 3, -2), (5, -2, -2) 
and (7, -7, -2) resonances, and they are connected in this high ec- 
centricity region, allowing bodies to switch from one resonance to 
another. 

The fact that some of the family members have stable orbits 
was the reason why we were able to apply backward integration 
method to estimate the age of the family. On the other hand, pres- 
ence of the chaos in the region occupied by the family, allows to 
use MCC in order to estimate the age. Using both methods, and 
combining the results, we found the age of Theobalda family to be 
6.9 ± 2.3 Myr. Given the very good agreement between results ob- 
tained with different methods as well as when applied to different 
groups, we believe this estimate is very robust. Thus, this is another 
family younger than 10 Myr. This result has several important im- 
plications, and some of them we mention bellow. 

The young asteroid familie s are also known to be source of 
solar system dust bands (see e.g. iGriin et al . 1985: Nesvornv et al.l 
[20031). The origin of three main dust bands is k nown, and they 
correspond to Karin, Veritas and Beagle family JNesvomv et al.l 
l2006bl : [Nesvornv et al"]|2008h . Also, the very young Emilkowal- 
ski fam ily is the most p robably source of incomplete dust band 
at 17° i JEspv et al.l [20091) . On the other hand, the origin of some 
less prominent bands, such as so called M/N dust band, is still 
not quite clear. The Theobalda family's young age, and its proper 
inclination of Ip ~ 14°. 3 suggest t hat it might be a possible 
source of M/N dust band {Ip ^ 15°) ( ISvkej[l990l) . On the other 
hand, this dust ba nd was linked to (170) Maria asteroid family 
I Reach et alj[ 19971). an d more recently to (1521) Seinajoki cluster 
I Nesvornv et alj[2003i) . In any case, dust band produced by such a 
young family, as Theobalda, should be observable. The size of its 
parent body also suggests that it should produce a prominent dust 
band. If this is not M/N dust band, then there must be another dust 
band which can be linked to this family. Alternatively, it should be 
explained why and how this dust band has disappeared. 

Theobalda asteroid family is located very close to the region 
where three (out of four) so -called main belt comets (MBCsr^l 
ha ve been disco vered (see e.g. [jewitt et alj2009f) . As was suggested 
bv [Hsiehl ( [20090 . it is possible that this kind of bodies can be found 
among the members of other young families, probably many of 
which waiting to be discovered. Being young and dominated by 



^^ The MBCs are bodi es with asteroid-lik e dynamical properties but comet 
like physical properties jHsieh et alj2004) . These are dynamically ordinary 
main-belt asteroids on which, probably, subsurface ice has recently been 
exposed e.g. because of a collision. 
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C-type asteroids, we believe Tlieobalda family is very good place 
to start. 
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